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The device of using "saturation arithmetic" to cope with adder over- 
flow in recursive digital filters has, for a number of years now, been 
known to yield stable operation when the filter is of second order and 
is linearly stable. Mitra has recently given examples to show that this 
happy situation does not prevail for higher order filters. Here we in- 
vestigate conditions on the filter coefficients which would guarantee 
stability for higher order filters using saturation arithmetic. We are 
only able to give sufficient conditions for stability. These conditions 
in their simplest form can be written as linear inequalities involving 
the coefficients of the filter. 



I. INTRODUCTION AND SUMMARY 

We shall be concerned with real nth order nonlinear difference 
equations of the form 

y(k + n) = f\zaiy(k + n-i)], k = 0, 1, 2, . . . (1) 

The variables y( • ) and the coefficients a, are real. The initial conditions 
y(j) t j = 0, 1, . . . , {n — 1) are arbitrary, subject only to the important 
condition \y(j) \ < 1. The function /(•) will be assumed here to have the 
form given in eq. (2): 

f(x)=x, \x\ < 1 

/(x) = sgnx, |x| >1 (2) 

This function models a method of handling overflow in the practical 
implementation of digital filters and in that literature is referred to as 
"saturation arithmetic." An important unsolved problem is the as- 
ymptotic stability of this undriven system. Specifically we would like 
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to describe the region in "a -space" (i.e. {a, j, i = 1 n) for which 

lim y n = 

for any initial conditions. We always assume that the a,-'s are already 
restricted so that linear stability holds. That is, if /(•) were replaced by 
the identity function the system would be stable. This is equivalent to 
all roots A, of the characteristic equation 

c(z)=z n - £ai*»-<-0 (3) 

satisfying | X,- 1 < 1.* Since {a,-} are real, the complex X, occur in conjugate 
pairs. 

The case n = 2 has special importance to a certain strategy for 
implementing digital filters and has been considered earlier. 1-3 It was 
shown for this case that eq. (1) is stable for saturation arithmetic 
whenever the system is linearly stable. Recently Mitra 4 has shown by 
example that a result of this generality does not hold for any n > 2. This 
surprising development has regenerated the author's interest in the 
problem in its own right. In addition, direct implementation of digital 
filters of the form [eq. (1)] for n > 2 is now of interest, so questions of 
stability must be answered. Saturation arithmetic seems to be an ex- 
perimentally favored procedure at the moment. 

Our main results, Theorems I and II, provide sufficient conditions that 
a given set of coefficients {a,-}? yield a stable filter with saturation 
arithmetic. Both theorems require one to test if a pair of linear ine- 
qualities in a set of variables can be satisfied when the latter lie in a hy- 
percube of (at most) dimension n. A finite algorithm which is sufficient 
to decide this question is given in Section IV. While we feel that Theorem 
II will give more powerful results (i.e., determine a larger stability region), 
Theorem I allows one to list a set of linear inequalities in the a,-, which, 
if any one is satisfied, would guarantee stability. This result is given as 
Corollary I. 

II. SOME LINEAR AND NONLINEAR THEORY 

If one were concerned with the linear version of eq. (1) [/(•) equal to 
the identity function] , the solutions would be 

y(k)=Zki\fk=0,l,... (4) 

;=i 



* Under this condition, stability of eq. (1) with |y,-| < 1 and 

LWi\<i 

is trivial, since the system is then linear. 
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where the A, are the characteristic roots of eq. (3), assumed to be distinct 
in writing the solution eq. (4). The constants /e, are determined by initial 
conditions, or equivalently, the values initially stored in the registers. 
Suppose we wish to solve for the fc, in terms of y (0), . . . , y(n - 1). This 
clearly requires the inversion of the matrix 



V(X) = 



1 
X? 



1 * 

An 
\} 

k n-l 



called a Vandermonde matrix. A brief discussion of these matrices is 
given in the Appendix, as well as some notation we shall use related to 
them. 

The linear difference equation can also be written in matrix form if 
we introduce the n -vector* 



Y(k) = 



y(k) 
y(k + 1) 



k = 0, 1, 2, . . . 



(5) 



y(k + n - 1) 

or, in component form y,(fc) = y(k + i - 1), i - 1, . . . , n. The "time" 
argument is indicated by the discrete index k. The equation-of-motion 



is then 



Y(k + 1) = AY(k) 



A = 






1 
















1 
































1 


a n 


Qn-1 




0-2 


Qi 



(6) 



(7) 



Expanding the determinant of the matrix A — XI (I being the identity 
matrix) by the last row, we obtain (except for a sign) the polynomial eq. 
(3), and, not surprisingly, the eigenvalues of A are the roots X, mentioned 
earlier. 

If these eigenvalues are distinct, a well-known theorem of algebra 
guarantees that there exists a nonsingular matrix P such that 



P~ l AP= A 
where A is simply the diagonal matrix of eigenvalues of A. 



(8) 



f Vectors and matrices will be denoted by capital letters. 
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Theorem: If the eigenvalues of A are distinct, then 

V~ 1 AV=A (9) 

where V is the Vandermonde matrix formed from the roots of the 

characteristic eq. (3). 

The choice P = V is not claimed to be unique. 

Proof: 



£ OijVjk = L aijH = 



-i _ I A * 



X' fc i<n-l 






The second line in the last member makes explicit use of the fact that 
\k is a root of the characteristic equation. Thus 

n -i / n \ n 

L vJi ( L QijVjk ) = E utt (XkVik) = ^k&ek 

j=i V=i / i=i 

as was to be shown. 

One reason for wishing to diagonalize A in the linear case is the simple 
form that the equation of motion takes. If we multiply eq. (6) by V -1 and 
perform the standard trick of inserting J = W" 1 after the A in eq. (6), 
we obtain 

Z(k + 1) = AZ(k) (10) 

where 

Z(k)=V- 1 Y(k) (11) 

Since A is diagonal, the solution for the ith component of Z is simply 

Zi (k) = \?Zi(0) (12) 

Turning now to nonlinear problems, we wish to summarize some re- 
sults from Liapunov stability theory, 7 without proofs, and without 
complete generality.* 

We are concerned with an autonomous (time independent) nonlinear 
difference equation 

Y(k + 1) = F[Y(k)] (13) 

where F is a nonlinear (or linear) vector function of the vector Y{k). 



* A. N. Willson was the first to explicitly apply Liapunov theory to the present problem 
for n = 2. 3 He has also attacked other stability questions for n = 2 with these methods in 
Ref. 8. 
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If 



lim Y(k) = 



then the system is called asymptotically stable. 

Theorem (Liapunov): If there exists a (strictly) positive definite qua- 
dratic form w[Y], such that, for any allowed n -vector Y, w[F(Y)] — w[ Y] 
< (strict inequality) then the system is asymptotically stable. + 

In other words, if we can find a positive definite quadratic form (of 
the state variables of the system) which is always decreasing as the mo- 
tion proceeds, then the motion must proceed to the origin. The function 
w[-\ is called a Liapunov function. 

In terms of /(•) and A the function F( ) is determined by 

[F(Y))i= [AY]i ;=l,...,(n-l) 
[AY] n , if \[AY] n \ <1 



[F(Y)) n = 



(14) 
sgn [A Y] n , otherwise 



As a simple example of a Liapunov function consider the linear case 
with nondegenerate eigenvalues, and again set Z = V -1 Y. Choose 

w=t\zi\ 2 (15) 

the zi being regarded as functions of the y (j '), j = 0, . . . , n — 1. We know 
that when Y — AY we have Z —» A.Z and so 

w-+il |X,| 2 |2,| 2 

1 = 1 

Since we assume | X, | < 1, strict decrease of w is assured. 

III. A SPECIAL LIAPUNOV FUNCTION 

For the nonlinear problem eq. (1), we shall, for a first pass, choose a 
w[-] whose form is inspired by the one just described. Noting from the 
Appendix 

n/-ivl (-D'" (/) (A) ^ , 1U (,) m 
[V M ^]/= 771 £ (-D J Pn-iih)yj 



u(\) e=o 

we shall single out the functional 



(_1 )i+n n U)(\) n-l 

X) e (-i)WMy(n-e-D (16) 



f For our problem any vector V is allowed that has components |y,-| < 1. We identify 
V with Y(0) and soy, = y{i - 1), i = 1 n. 
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xi = E 1 (-l)<pf(\)y n -e = E 1 (-D'plP<A)y<n -£-l) (17) 

for special attention. Clearly the equation of motion for the Z variables 
(10) implies that under the substitution Y -*■ A Y we have * t -»- X,-*,-. We 
choose 

w[Y)=t\xi\ 2 (18) 

i=i 

where the x u via eq. (17), are regarded as functions of y;, i = 1, . . . , n, 
the components of V. Of course we always have |y, | < 1. 

In order to investigate the consequences of the (sufficient) stability 
condition 

w[F(Y)]-w[Y]<0 (19) 

we note that under Y -» AY we have x, — x\ L) (L stands for linear) 
where 

x\ L) = "e (-lKp]°(A)y(n - *) (20) 

and 

y(n) = E a^y(n - *) (21) 

^=i 

Finally the nonlinear (NL) "version" of eq. (20) is 



^J«f»*M»>l«i (22) 

sgny(n) + E (-lKp^(X)y(n - *) if |y(n)| > 1 

i 

where we made use of the definition of F(-), and the fact that po (£) { • ) 
= 1. Then 



w[F(Y)]=Z \xl NL) \ 2 (23) 



i=i 



We have already noted that if \y{n)\ < 1 we have a linear iteration 
[F(Y) = AY] and 

w[F(Y)] - w[Y] = w[AY] - w[Y] 



-E |x.-| 2 l*.-| 2 -E l*.-| 2 <o (24) 

i'=l i-l 

Hence we will only be concerned with :c j NL) when \y(n) | > 1. In this case, 
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w[F{Y)] is a function of y(l),y(2), . . . ,y(n - 1), while the condition y(n) 
> 1 involves y(0) as well.* 

We shall not use the stability condition in the form of eq. (19), but 
instead note that since 



El*Pl a -El7i| 2 W a <i:W* (25) 

1=1 



the condition 



E |*! NL) | 2 - £ |x/ L) | 2 <0wheny(n)>l (26) 

i = \ j' = l 

is sufficient for stability. At the price of losing some power in the method 
we shall see momentarily that we have gained considerably in analytic 
simplicity. For convenience define 

d = - "e (-l)'y(n - £W{\) (27) 

£=l 

so that [when y(n) > 1] 

x\ L) = y(n)-ci 
X (NU = i _ c . (28) 

Then 

E |*P| 2 - E k^ L | 2 = E (y(n) - l)(y(n) + i- Ci - c*) 
i i i=i 

= (y(n) - 1) \n(y(n) + 1) - 2 £ a j (29) 
We have used the fact that complex X, occur in conjugate pairs and so 

i i 

Thus if the inequalities 

2 tci-n(l + y(n))>0 
i=i 

y(n) > 1 (30) 

have no solution in the rc-cube |y (0 1 < 1, i = 0, 1, . . . , n — 1 the nonlinear 
equation (1) is stable. More explicitly by using the definition of the c, 
[eq. (27)] Lemma I in the Appendix coupled with (67) to express 2f =1 c, 
in terms of the a,, and using finally eq. (21), we have: 



* By symmetry of the problem, |y (n) | > 1 is here and henceforth replaced with y (n) >1. 
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Theorem I: If the two inequalities 

Z (n - 2£)a e y(n -£)>n (31) 

£=l 

Z a £ y(n - £) > 1, (32) 

cannot be simultaneously satisfied by some set of y (i) in the cube \y{i)\ 
< 1, i = 0, 1, . . . , n — 1, then the system (1) is stable for saturation 
arithmetic. 

We note that if the inequalities are satisfied, no conclusion is 
drawn. 

A systematic algorithm for checking the above inequalities is given 
in the next section. Here we deduce the simple: 

Corollary I: If the coefficients ae satisfy 

Z \a £ \\n-2£\ <n (33) 

then eq. (1) is stable. If the coefficients a\ satisfy 

E \a e \\k-e\ <k (34) 

e=\ 

for at least one k, [n/2] < k < n, eq. (1) is stable.* 

The first inequality follows immediately from the first inequality in 
Theorem I plus the fact that \y(i) | < 1, all i. The remaining inequalities 
follow from the observation that in eq. (31), the coefficients of a^ have 
the same sign for £ < [n/2], whereas for £ > [n/2] they have opposite 
signs. Thus if eqs. (31-32) have a simultaneous solution, so do 

Z \a £ \(n-2£)+ Z (n - 2£)a £ y(n - £) > n (35) 

l<£<[n/2] £>[n/2] 

Z \at\+ Z a e y(n-£)>l (36) 

£<[n/2] £>[n/2] 

where the above is obtained by setting y(n — £) = sgn a £ , 1 < £ < [n/2]. 
If, for k > [n/2], we multiply the second inequality by (2k — n) > and 
add the result to the first we obtain 

Z \a e \(k-£)+ Z (k-£)a £ y(n-£)>k (37) 

£<[n/2] £>[n/2] 

which, if \y(i)\ < 1, cannot possibly be satisfied when 



* The notation [x] denotes the integer part of x. 
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£ |a,||*-*| <k (38) 

There is one important point to be noted here. Our results have only 
been proven for the case of nondegenerate eigenvalues. For degenerate 
eigenvalues, the Liapunov function that we have chosen in this section 
is not strictly positive definite. We have in fact constructed Liapunov 
functions which are specifically designed to handle the degenerate case. 
Using them, we have proven that the conclusions are still true for de- 
generate eigenvalues. We believe that the form of the results, being 
simple conditions on the a\s, will allow the reader to readily accept that 
they are true in general. Since our proof of the extension is long and out 
of proportion to its importance, we have chosen to omit it. 

IV. HYPERPLANE ALGORITHM 

Let (6/|f=i and {c,-}?=i, £ and 77 denote fixed constants. We wish to de- 
termine when it is possible to simultaneously satisfy the inequalities 

E Zibi > £ 
k 

1=1 
\zi\ <1 1 = 1,2 k (39) 

The dimensionality of the problem is immediately reduced if sgn bj = 
sgn Cj for some j since we may immediately take zj = sgn bj. It is im- 
portant to note that we assume this to have been done and therefore 
assume 6,c, < 0, 1 < i < k. 

Lemma: If the simultaneous inequalities eq. (39) are satisfied then there 
exists zi, and ay, 1 < j < k, such that 

L 2^ > ^ 

1=1 

k 

L ZiCi > rj 
1 = 1 

\tj\ <l,\2t\ =la\\£*j (40) 

In other words all but possibly one of the coordinates may be given values 
±1. 

This is geometrically evident if/? = 2. If k > 2 one need only consider 
the Zi variables two at a time, always applying the Lemma for k = 2. 
Eventually all but perhaps one of the 2, will have value ±1. 

Continuing with the description of the algorithm, let E be any k -vector 
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having components t£ = ±1; there are 2 k such E vectors. Choose a/ and 
let's test if indeed it is the j of the Lemma. Let z~e = te,£ ?* j. Then if eq. 
(40) has this solution we have 

g - 2'b £ e £ ^ . ^ i? - 2'e £ c £ 

< z; < if bj > (41a) 

bj Cj 

or 

£^^<2.<i^f**<if 6 ,<0 (41b) 

If these inequalities are consistent (i.e., the upper bound is at least as 
big as the lower bound) and if they can be satisfied by some Zj, \zj\ < 
1 we are done — the simultaneous inequalities are satisfied. If not, try 
another E vector, or another;. For a given; there are 2 k ~ l E vectors to 
try. Hence after at most k-2 k ~ l such attempts we have exhausted all 
things that need to be checked, and checking the inequalities is, it has 
turned out, a finite procedure. 

This procedure is not only applicable to Theorem I, but also to The- 
orem II occurring in Section V. 



V. ANOTHER LIAPUNOV FUNCTION 

We have already noted that the entire sequence \yi\o is determined 
by the first n elements. For our second choice of the Liapunov function 
we chose the expression for the energy in the remainder of the sequence 
for the linear problem: 

w[Y] = flyk (linear case) (42) 

n 

The right member is regarded as a positive definite quadratic form in 
Y(0). In the linear case we also have, numerically, 

w[AY] = ±yl (43) 

n+i 

which is smaller than w[Y] by y 2 (n). Thus this w[-] doesn't necessarily 
decrease after every iteration and thus it is not strictly a Liapunov 
function. However after at most n iterations it must decrease (unless all 
yi = 0) and the effect will be the same. We shall have stability if we can 
show whenever y (n) > 1, that 

w[F(Y)]-w[Y]<0 (44) 

or, equivalently 

w[F(Y)]-w[AY]<yl (45) 
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We begin by writing down the generating function for the sequence 
y(n), y(n + 1), . . . , when linear theory holds. By definition 

H(z) = £ y(n + k)z k (46) 

Using standard linear techniques we calculate from eq. (1) and its initial 
conditions 

n n—j 

E y(n-j) E aj+sZ s 

H(z)=^ — ^ (47) 

- E a s z° 

s=0 

where we have arbitrarily defined ao = —1. We note the characteristic 
polynomial 

c(z) - - E aiz *"* 

has the same modulus as the denominator ofH(z) when \z\ = 1 (since 
the a, are real). 
We next note that 

£yl = ^- r\H(z = e")\*d6 (48) 

To express this as a quadratic form introduce, for 0< |s — t\ < (n — 1) 
the integrals (z = e ie ) 

1 ~« z s-t dd = J_ pe ep, (, - t )6 de m 

st 27T Jn I „ 2 2tt J-.i « 2 

E a s z s E a„z s 

ls=0 I s=0 

Also, as a contour integral around the unit circle, we have 



s ' 2iriJP (2a s 2 s )(2a t 2"-0 2 



(50) 



The first form of the integrals shows they are real and I st = I ts - Also in- 
troduce the real symmetric matrix 

n—j n—k 

Hjk = E E aj+sOk+tlst, J, k = 1, . . . , n (51) 

s=0 * = 

Note for j = k, Hjj > since the right side of eq. (51) is then a Teoplitz 
form with positive spectral function. Then with these notations simply 
using eq. (47) to expand the integral in eq. (48) yields 

w[Y] = £ yl = t y(n -j)y(n - k)H jk (52) 

k=n j,k = l 
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or 

w[AY) = £ y\ - £ y(n + 1 - ;)y(n + 1 - fe)Jfy (53) 

When y n > 1, h>[F( Y)] has the same form as eq. (53) except thaty n is 
replaced by unity. Thus in evaluating u;[F( Y)] — w[A Y] [by using eq. 
(53)] most of the quadratic terms will cancel. Doing this and a few minor 
manipulations, the criterion for stability will read that we must have 

\y(n) + l)H n + 2 £ H ljy (n + 1 - j) > ^L (54) 

y=2 y(n) - 1 

whenever 

yiXt)- T.a jy (n-j)>l (55) 

If we define Hi n +\ = 0, we may write this as: 

Theorem II: If there are no simultaneous solutions to the inequalities 

£ l-a&u ~ 2ffi,»+i]y(n " £ Hn + f^ : (56) 

i-i y(n) - 1 

y(n)= £ ajy(n-j)>l (57) 

7 = 1 

|y(i)|<l i = 0,l,...,n-l (58) 

then the difference eq. (1) is stable with saturation arithmetic. 

To convert this into a hyperplane problem (discussed in Section IV) 
the nonlinear term y 2 (n)/(y(n) — 1) may be dropped or replaced by the 
value four (since 

x 2 

>4 

x -1 

when x > 1). If we drop this term, the resulting stability has the physical 
interpretation that at any time the energy in the remaining tail of the 
nonlinear response is less than or equal to the corresponding energy for 
the linear problem, regardless of the previous state. That is, if measured 
by energy, the nonlinear undriven response dies off at least as fast as the 
linear response for any initial conditions. 

We also note that whenever the Liapunov function w = Y^HY [eq. 
(52)] drops to the value unity the system behaves as a linear one from 
there on (no future yk will exceed unity). Several bounds for this quantity 
may be given. Since 

2 



E a .s- 2 ' s 

5-0 



= n |i-*x £ |2>n(i-|A f |) 2 
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we have 



w 



-— C" dd\H(z)\ 2 < max \H(z)\ 

2,1T %) — -k z=e' n 



= max 





-/ 


1-J 

) E 

s=0 


«/+■* s 




n 

E 




a s 2 s 





n-j 



E |y(« -;)| E !««/+« I 

>=1 s-0 



E |o;| 

y=i 



n E|y(0| 

j=0 



(59) 



n (i- |x ( |>2 
i=i 



n (i - ix ( |) 2 

j=i 



VI. EXAMPLES 

The few simple examples of this section will shed light on the two 
methods we have given. When n = 2 the two inequalities of Corollary 
I are simply | a 2 | < 1, |ai| < 2. Since linear stability implies |ai| = |Xi 
+ X2I < 2, |a2| = l^iM < 1 1 we see that torn = 2 linear stability implies 
stability with saturation arithmetic. We have already mentioned this 
is not true for n > 2. Mitra has constructed a counter example using the 
degenerate case Xi = X2 = . . . = X n = 7. For n = 3 he finds oscillations 
if 1 7 1 > 0.858, although if 1 7 1 is smaller than this stability is not implied. 
If we consider the second inequality of Corollary I forn = 3, k = 2, a.\ = 
37, ct2 = — 37 2 , 03 = 7 3 we have stability if 

3|7l + M 3 <2 

or I7I < 0.596. No better result is obtained for this case by a complete 
use of Theorem I. 

On the other hand if we apply the criterion of Theorem II (neglecting 
the nonlinear term) with the algorithm of Section IV, we find stability 
if I7I < 0.71. Insignificant improvement would be obtained here if we 
had also included the nonlinear term. The application of Theorem II to 
the present case was sufficiently simple so that the calculation could be 
done by hand. The integrals were done exactly to give 



H u = 



X 1 



(1 - \2)5 



[9 - 9X 2 + 10X 4 - 5X 6 + X 8 ] 



^12= 3A L (3 + X 2 ) 



Hia — 



(1-X 2 ) 5 
3X 4 

(1 - X2)5 



(1 + X 2 ) 



(60) 
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Clearly in general the integrals would have to be done numerically. The 
fact that the nonlinear term did not contribute significantly is due to the 
combination of facts that at the critical value for X,y(n = 3) is not ex- 
tremely close to one and also Hu is large, due to its denominator. 
It should be pointed out that 

n 

can be expressed using the solution eq. (4) for the linear problem as 

Y,yl = Y+iV-^rv-tY (61) 

n 

where 

Tij= , .' (62) 

However this explicit form is only true for the nondegenerate case and, 
although all limits exist as Xi -*■ X2, etc., the expression would probably 
not be suitable for numerical computation when eigenvalues are close 
to being degenerate. 

Another example of limiting misbehavior being only apparent is that 
in a similar manner one could compute that 

H(z)= E- J - r (V- 1 Y) i (63) 

;=i 1 - zX," 

Individual terms in this expression are badly behaved as, for example, 
if Xi « X2, but the alternate form eq. (47) shows everything is well be- 
haved in the limit. 

If we return to Corollary I applied to n = 3 when (Xi, X2, X3) = p(i, —i, 
1), < p < 1, we see that the inequality 

|ai| + |a 3 | <2 

is sufficient to guarantee that for filter poles in these relative position 
saturation arithmetic will give a stable filter for any p, all the way out 
to the boundary of linear stability. 

Based on these examples we feel that Theorem II is the more powerful 
method although the simpler Corollary I can yield considerable infor- 
mation for particular cases. 

Finally, we note that an important investigation on the present 
problem has just been completed by Mitra, 9 resulting in different sta- 
bility criteria from those presented here. Mitra's results will give a 
polynomial type criterion for absence of periodic oscillations. These 
results in themselves do not prove stability in that Ref. 9 does not pre- 
clude unending periodic outputs with no input. However, we take the 
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liberty of mentioning that Mitra has extended the proof to include sta- 
bility and thus the results of Ref. 9 may be taken as proving the same 
type of stability as discussed here. Another comparison with Ref. 9 in- 
volves the size of the stability region in "tap-space" which the two 
methods give. Neither Mitra's criterion nor ours can claim to describe 
the largest stability region. Also it does not even seem possible at this 
stage to give theoretical arguments to decide if one of the methods is 
always superior in this respect. However, several examples indicate that 
the region determined by Mitra's criterion is larger. Assuming this to 
be the case in general, an effective practical procedure would be to first 
test for stability using our simple Corollary I, and if this fails, apply 
Mitra's polynomial test. 



APPENDIX 

The Vandermonde Matrix and Symmetric Polynomials 

The a denote an ordered set of n complex members a,-, i.e., a\, 012, . 
a n . By the Vandermonde matrix V(a), we shall mean the matrix 



V(a) = 



«i 

4 



«r 



« 2 
Oil 



or 1 



as 

4 



ar 1 



otn 

„.2 



..n-1 



(64) 



This can be written [ V(a)],y = a) l ,i,j— I, .... n. We let v(a) = det 
V(a), and it is known 5 that 



v(a) = U{aj - or,-) 
where the product extends over all i, j satisfying 

1 < i < j < n 



(65) 



If a, 5^ aj for i 5^ j then the inverse of V(a) exists, and is known. Before 
giving its structure, we wish to list some facts concerning some special 
symmetric polynomials. 6 

Definition: The ^th elementary symmetric function of n -variables (£ 
— 1, 2, . . . , n) is the sum of all formally distinct products of the variables 
taken £ at a time. We also define po - 1. 

For example if n = 3 we have 
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Pote) = 1 

Pi(a) = ai + a 2 + oc 3 

p 2 (a) = ai<x 2 + «ia 3 + a 2 «3 

P3(a) = aia 2 a 3 (66) 

A well known theorem of algebra states that any symmetric polynomial 
in all the a's can be uniquely written as a polynomial in the quantities 
Pi(a),i = 0,...,n. 

Note that in the characteristic polynomial eq. (3) we have 

ai=(-iy +l Pi (\) (67) 

where the X,-, i = 1, . . . , n, are the roots of eq. (3). Thus the theorem just 
stated says that any symmetric polynomial in the roots of a polynomial 
can be expressed as a polynomial in the coefficients of the equation 
(rather than a complicated function as would be required to express an 
individual root). 

We shall use the notation p^ l) (a), £ = 0, . . . , n — 1, to denote the ^th 
elementary symmetric function formed from the (n — 1) ordered vari- 
ables ai,..., otj-i, oti+i, . . . , a n . Likewise d ( 'H«) denotes the determi- 
nant of the corresponding (n — 1) X (n — 1) Vandermonde matrix. 

Theorem: 

lV-H a )] U . <- 1)i+ X Ka) Piik(°) (68) 

u(a) 

i, j = 1, . . . , n 

Proof: We have 

Qkj= L, [V{a)]ki[V Ha)\ij- £ a* — l (69) 

From eq. (65) and the definition of u (l) (a) we see that 
u <«>(«) 1 



v(a) (-i)n-ij-j ( a . _ a£ ) 



(70) 



and hence eq. (69) is 

^,.~ af->g>, ■(«)(-!)-/ (?1) 

i=i n («i _ a*) 
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and note that* 



Thus 



q(x) = £ qkjXJ 



£ pHljMi-l^-JxJ = L (x - at) (72) 

y'=l ^j" 



IT (* - at) 

£ QkjX'- 1 = L af" 1 Tf- " (73) 

;-i '=1 n («i — «*) 

has value a£, -1 when x = a m , m = 1, . . . , n. From this it follows that q(x) 
= x* -1 , so that qkj = &kj> which we were to prove. 

We leave it to the reader to convince himself of the following: 

Lemma I: 

t P { n%M = jp n -j(a) ;>0 (74) 

e=\ 

Lemma II: 

p ; (X)-pj' ) (X) = X t pj'2 1 (X) j = l,2,...,n-l (75) 
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* The right member of (72) is a polynomial expressed directly in terms of roots, the left 
member, via (67), the polynomial expressed directly in terms of coefficients. 
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